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We present a method to detect the unstable periodic orbits of a multidimensional 
chaotic dynamical system. Our approach allows us to locate in an efficient way the 
unstable cycles of, in principle, arbitrary length with a high accuracy. Based on a 
universal set of linear transformations the originally unstable periodic orbits are trans- 
formed into stable ones and can consequently be detected and analyzed easily. This 
^ ' method is applicable to dynamical systems of any dimension and requires no preknowl- 

edges with respect to the solutions of the original chaotic system. As an example of 
application of our method we investigate the Ikeda attractor in some detail. 
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I ■ I. INTRODUCTION 

o 

The finding that the unstable periodic orbits provide a skeleton for the organization of the very 
complex chaotic dynamics can be considered as one of the major advances in our understanding of 
the behaviour of nonlinear dynamical systems during the past ten years. Many physical quantities of 
these systems can be determined by knowing the positions/properties of the unstable periodic orbits 
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living in the chaotic sea. For strange attractors belonging to dissipative systems the quantitative 
description of their structural properties in terms of unstable periodic orbits has been successfully 
used for several dynamical systems This reflects the importance of the cycles for the analysis 

and decoding of the dynamics on the attractor. On the other hand knowing the positions of the cycles 
in phase space one can use them to control chaotic dynamical systems Moreover the quantum 
mechanical properties of classically chaotic conservative systems possess in the semiclassical regime a 
series expansion with respect to the lengths and the stability coefficients of the periodic orbits (see 
ref. |t]] and references therein). Since chaotic behaviour is an intrinsic property of many dynamical 
systems periodic orbit theory possesses numerous applications in different areas of physics. 

Assuming that the cycles represent the skeleton for the chaotic dynamical evolution of dissipative 
systems one can use an appropriate series expansion in terms of their period to determine observables 
like fractal dimensions, average Lyapunov exponents, entropies or the invariant measure of the corre- 
sponding strange attractor . The convergence properties of this expansion allow us to test the 
applicability of current periodic orbit theory for the properties of dissipative dynamical systems. The 
most important and at the same time most difficult step of this procedure is to find the location of the 
unstable cycles in the chaotic sea: in spite of the fact that their number is growing exponentially with 
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increasing period and that they are dense in the chaotic sea the unstable cycles are of measure zero 
whereas the chaotic orbits are (for a fully chaotic system) of measure one. However, once the location 
of the orbits is known it is straightforward to extract their properties from the underlying equations 
of motion. It is therefore not surprising that a number of recent works deal with the development of 
efficient methods and strategies for the detection of the periodic orbits |]lO|"|l^|. 

The application of the Newton-Raphson method to find the unstable cycles requires a good guess 
for the starting point of the iterative procedure. For cycles of higher period this is only possible if one 
uses a very fine grid of initial conditions on the attractor. Consequently this approach involves a huge 
numerical effort for more than two dimensions and is therefore in general not feasible. Alternative 
methods can be found in the literature but their applicability is limited to special low-dimensional 
systems @|[|]l§. 

The purpose of the present investigation is to establish a generally applicable as well as reliable and 
accurate method in order to determine the unstable periodic orbits in a chaotic dynamical system 
H]. The basic idea of our approach is to use an universal set of linear transformations in order to 
transform the unstable periodic orbits of the original system into stable periodic orbits which occupy 
the same positions in coordinate space. The stable periodic orbits of the transformed systems can 
then be found by simply iterating these dissipative transformed systems. It can be shown that such 
transformations exist in general and possess very restrictive and simple forms as well as an appealing 
geometrical interpretation in terms of a vector field which is organized through the positions of the 
cycles. To demonstrate the efficiency of our method we apply it to the attractor of the two dimensional 
Ikeda map |l5| which describes the dynamics of an optical field in a ring cavity. We locate the unstable 
periodic orbits up to period 13 of the Ikeda strange attractor. Subsequently we derive their stability 
coefficients and use the expansions according to periodic orbit theory in order to determine the average 
Lyapunov exponents, the fractal dimension as well as the topological entropy of the Ikeda attractor. 

The paper is organized as follows : In section 2 we present the general theoretical framework 
of our method and discuss in more detail the case of a two dimensional system. Futhermore we 
present a universal algorithm for determining the unstable periodic orbits of a multidimensional chaotic 
dynamical system. In section 3 we use this algorithm to locate the unstable periodic orbits of the Ikeda 
attractor. From their location we derive the corresponding stability coefficients which can be used to 
determine the average Lyapunov exponents, the fractal dimension and the topological entropy of the 
attractor. In section 4 we summarize the main aspects of our approach and briefly report on possible 
applications as well as perspectives related to the general problem of solving nonlinear equations. 



II. THE THEORETICAL FRAMEWORK 

Let us consider a TV-dimensional discrete fully chaotic dynamical system given by: 

u : n+i = M) (i) 

U, being fully chaotic, possesses only unstable fixed points (FP). Since points of a periodic orbit of 
period p are FP of the p — th iterate /' p ' the term FP stands in the following for periodic orbits 
in general, i.e. for orbits of any period. We thereby have to replace / in eq.(I) by f^ p \ Our goal 
is to construct from the map U in eq.(|l|) different dynamical systems {Sk\k = 1, ...,M} (see below) 
which possess FP at the same positions as U but instead of being unstable they have become stable 
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in the dynamical systems {<Sfc}. The corresponding transformations {Lk : U — > Sk} are required to 
preserve the number of FP, i.e. no additional orbits should be created. The transformations {Lk} 
therefore change the stability properties but not the location of the FP. The index k stands for the 
possibility that different orbits may require different transformations for their stabilization. However, 
as we shall see below, each type of transformation will stabilize not only a single FP but a whole 
rather general class of infinitely many FP ranging up to arbitrarily high periods. If we succeed with 
our plan then the search for the positions of the FP of the system U becomes a manageable task: 
because of the stability of the FP in the transformed dissipative systems {Sk} each trajectory of Sk 
after some iterations closely approaches a FP fp ■ Per construction fp is then also a FP of the system 
U and we therefore know the position of the FP in the original system U. 

To fulfill the requirement of the one to one correspondence between the FP of U and Sk the 
transformation Lk should in general be linear. Consequently Sk takes on the following appearance: 

S k ■■ n+i = n + Afc(/(fi) - n) (2) 

where is a constant invertible N x N matrix. The definition (Q) satisfies the one to one corre- 
spondence of the FP of U and those of Sk- If f = rp is a FP of U then the parenthesis on the right 
hand side of eq.(Q) vanish and therefore rp is also a FP of Sk- On the other hand if fp is a FP of Sk 
and since is nonsingular the parenthesis on the right hand side of eq.{^) must be equal to zero for 
f = rp which implies that rp is also a FP of U. Thus the dynamical laws U and Sk possess FP at 
identical positions in space. 

Our next step is to stabilize the FP of the transformed systems Sk by suitable choices for A*,. 
Different unstable FP of the map U arc then, in general, stable in different transformed dynamical 
laws Sk- However, as we shall see below, the number of matrices {A&} necessary to achieve stabilization 
of all FP is very small. It turns out that if the absolute values of the elements of the matrices Aj, 
are sufficiently small (|Ay | <C 1, i,j — 1, »,N) then there exists a universal set of very restrictive 
matrices such that at least one matrix A& belonging to this set transforms (via eq.(j^)) a given unstable 
FP of U to a stable FP of the corresponding map Sk ■ In order to determine this set of matrices let 
us consider the stability matrices Tjj and T$ k of U and Sk which obey the following relation: 

M k : T Sk = l + A fc (2V-l) (3) 

For Tjj we assume that it is real, invertible and diagonalizable. Since rp is an unstable FP at least one 
of the eigenvalues of Txj at rp must possess an absolute value greater than one. In order to stabilize 
fp we proceed in two steps: first we use the parametrization (Afc)y = [XCkjij with I ^> A > 
and Cij = 0(1). The matrix Cy has to be chosen such that the real parts of all eigenvalues of the 
matrix Ck • (Tjj — 1) are negative. If this is achieved then the next step is to use a sufficiently small 
value for the parameter A, such that the eigenvalues of the matrix Tg h = 1 + XCk{Tjj — 1) have 
absolute values less than one. It can be shown that this is always possible if A is sufficiently small: 
for sufficiently small values of A the leading contribution with respect to A in the absolute values of 
the, in general complex, eigenvalues is of order A 1 and emerges from the real part of the eigenvalues 
of the transformed system. Consequently we have transformed the unstable FP fp of U via Lk into 
a stable FP of Sk thereby keeping its position fixed. 

Next let us derive possible sets of matrices {Ck} which, according to the above discussion, transform 
the signs of the real part of the eigenvalues of Tjj — 1 and consequently allow us to stabilize any 
configuration of unstable FP. First we provide a rather general set of matrices and subsequently 
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we will specialize to a much more restrictive and simpler set. For any given Tjj it is possible to 
find an involutory matrix Ck defined through C\ = 1 such that A — Ck{Tu — 1) has eigenvalues 
with negative real parts. The proof reads as follows: Assuming the diagonalizability of the matrix 
B = Tjj — 1 there exists a similarity transformation P such that Bn = P~ X BP is diagonal. We 
have then: 



We can therefore always choose a diagonal matrix Cnk = P~ x CkP with elements ±1 on its diagonal 
such that Cnk ■ Bn has eigenvalues with negative real part. Such a matrix Cnk is according to its 
definition involutory. Due to the invariance of the eigenvalue spectrum of a matrix with respect to 
similarity transformations, the matrix A has then also eigenvalues with negative real part. In addition 
if Cnk is involutory then any matrix resulting from Cnk through similarity transformations, and in 
particular Ck, is also involutory. We can therefore always find an involutory matrix Ck given by 
Ck = PCnkP 1 such that Ck{Tu — 1) possesses eigenvalues with negative real parts. 

The set of involutory matrices is for our purposes rather general. For reasons of applicability as well 
as practical efficiency of our method it is therefore desirable to find a smaller set of transformation 
matrices. In fact it turns out that a much more restricted form for the matrices Ck is sufficient 
in order to achieve stabilization of any FP via the transformations {Lk}, namely all the matrices 
corresponding to special reflections and rotations in space. The elements of this new class of matrices 
{Ck} are C\V € {0, ±1} and each row or column contains only one element which is different from 
zero. The matrices Ck are therefore orthogonal. The number a at of such matrices in N dimensional 
space is given by cin = Nl 2 N . 

To illuminate our approach we discuss in the following the case N = 2 in more detail. In this case 
U is given by: 



and represents a fully chaotic 2-D map with a 2 x 2 stability matrix Tjj. Let us denote with pi.2 the 
eigenvalues of Tjj. Being fully chaotic U possesses only hyperbolic FP. The choice of the matrix Ck 
appropriate for the stabilization of a particular FP depends on the class to which the hyperbolic FP 
belongs. According to the above discussion there are a-i = 8 possible matrices {Ck} as candidates to 
be used for this stabilization process. With some algebra one can find the matrices which achieve the 
stabilization of the different types of unstable FP: 

• For hyperbolic FP with reflection and p\ < — 1, — 1</52<0=>Ci = 1 

• For hyperbolic FP with reflection and pi < — 1, 0<p2<l^-C'i = l 

• For hyperbolic FP with reflection and p\ > 1, — 1 < p% < or hyperbolic FP without reflection 
(pi > 1, < p2 < 1) we have to distinguish between the following three cases: 



P 1 AP = P 1 C k P B D 



Xi+i = f{xi,yi) 

Vi+i = g{xi,Vi) 



(4) 
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With the above choices for Ck and sufficiently small values of the parameter A we can stabilize any 
hyperbolic FP of a chaotic 2-D system. The only exception where our approach does not work is the 
case of a parabolic FP, i.e. for parallel flow (pi = 1 V pi = 1). However this case is of no interest to 
the present investigation. 

It is a major advantage of the present approach that a single matrix {C^} is responsible for the 

stabilization of an infinite number of unstable FP belonging to periodic orbits of arbitrarily high 

periods. For example: the matrix C\ stabilizes hyperbolic FP (periodic orbits) with reflection with 

respect to both invariant manifolds, independent of their period. We remark that only 5 of the 

8 matrices Ck are necessary to stabilize all hyperbolic FP in two dimensions. The matrices C4 

and Cg typically stabilize no or only a very few additional FP which is due to the fact that the 
df dg 

condition — — \pp = — \pp is met only very rarely. A similar statement holds also for higher 
ox dy 

dimensional systems: only a subset of the above given class of ajv matrices is actually necessary in 
order to achieve the stabilization of any given FP. According to the above we have shown a one way 
criterium with respect to the stabilization of the FP. For example: any hyperbolic FP with reflection 
(pi < — 1, — 1 < P2 < 0) in two dimensions can in particular be stabilized by using the matrix C\ but 
might also be stabilized by using other matrices of the set C/. . The explicit statements saying which 
matrices stabilize which types of FP become increasingly more complex with increasing dimension of 
the dynamical system under investigation. 

Having presented the basic features of our method we discuss now briefly how it can be used to 
detect the unstable FP of a given TV-dimensional discrete fully chaotic dynamical system. Using eq. (||) 
we transform the given dynamical system U into a new system Sk- For the matrix = AC& we use 
a sufficiently small value of A and C\ from {Ck}, k = 1, ..,azy Iterating an arbitrary initial point 
Fq with Si from eq.(^|) forward in time we observe the following characteristic behaviour: either the 
trajectory runs to a stabilized FP fp with steps of continuously decreasing size or, if stabilization is 
not achieved using Ci, the trajectory chaotically evolves on the attractor or escapes to infinity. In 
order to get all the FP of U this procedure has to be repeated for a representative set of starting 
points which covers in a crude way (see below) the phase space of the system. Subsequently we 
perform the same procedure for the next matrix C2 of the set Ck, etc... until all the matrices {Ck} 
have been used. Following this procedure we obtain the complete set of FP of U. In order to find 
the unstable FP of period i we simply have to replace / in eq.(|J) by its i-th iterate /\ One must 
however be careful with respect to the choice of the value of A in A& because it is indirectly related 
to the stability coefficients of the desired cycle. We therefore have to use decreasing values of A with 
increasing periods of the cycles we would like to detect. However later on in the present section we 
will provide a continuous formulation of the transformed systems which is independent of the value 
of the parameter A. Due to the fact that different kinds of FP (e.g. hyperbolic FP with and different 
kinds without reflection) are stabilized by different matrices Ck our stabilization procedure offers the 
possibility to distinguish between the different types of FP. 

An important advantage of our stabilization method is its global character. Even points lying far 
from the linear neighbourhood of the stabilized FP are attracted to it after a few iterations of the 
transformed dynamical law. To each FP of a certain period we can assign its basin of attraction 
consisting of the set of starting points which converge towards the stabilized FP if we iterate them 
with the transformed dynamical system Sk- With increasing period the typical volume of a basin 
of attraction becomes increasingly smaller and the basins form a geometrically very complex and 
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interwoven network which covers the attractor. 

Using the above approach to detect the unstable periodic orbits of a given chaotic dynamical system 
it turns out that the number of starting points needed to obtain the periodic orbits of a given period 
on the (closure of the) attractor is only a few times more than the number of cycle points themselves. 
Suitable starting points for the application of the transformed dynamical law can be obtained by 
using, for example, a chaotic trajectory on the attractor itself or more refined methods which cover 
the attractor in a more systematic manner. We used a simple empirical strategy in order to seek for 
the periodic orbits of a given period. Assuming that a certain number of periodic orbits of period p 
have been found we iterate a multiple of the starting points used in the previous run. If no additional 
periodic orbits show up we take the set of detected orbits as complete. We emphasize that this is by 
no means a proof for the completeness of the set of detected periodic orbits. However it turns out that 
for dynamical systems for which the periodic orbits are known, like for example the Henon map, the 
above procedure yields all periodic orbits. Because of the fact that each FP possesses a finite basin 
of attraction (see below and Fig. 2) the above method turns out to be rather insensitive with respect 
to the specific choice of the set of starting points covering the attractor in a crude and large-meshed 
manner. 

The transformed dynamical law (||) possesses an appealing geometrical interpretation. Let {fj , j — 
1, ..,p} be a trajectory of the original system U. At each point of the trajectory we define a vector field 
Vuifj) = fj + i — fj. The corresponding transformation Lk represents then a special reflection/rotation 
of each vector Vo (fj ) combined with a subsequent scale transformation of its length with the factor 
A. Using eqs.(l,2) we obtain for the vector field Vs k of the transformed map 

V Sk (f) = A k Vu(r) (5) 

The main feature of the new vector field Vs k is its global organization around those FP which have 
been stabilized. The flow of the vector field Vs k is centered and organized around the positions of the 
stable FP which represent sinks/sources of this vector field. To illustrate these properties we show in 
Fig.l selected vector fields belonging to a chaotic trajectory on the attractor of the Ikeda map (see also 
section 3). In Fig. 1(a) the vector field Vs k for the first iterate f 1 is illustrated. We use a trajectory 
of 200 points and the stabilization matrix Afc = 0.1 • C±. The global organization of the flow towards 
the FP fp — (0.5328,0.2469) whose position is indicated in Fig.l by a cross is evident. Obviously 
this property is not restricted or specific for the linear neighbourhood of the FP but represents a 
global feature of the dynamical system. Starting with any point on the attractor the trajectory of the 
transformed map moves immediately towards the FP and yields with increasing number of iterations 
increasingly accurate values for its position. In Fig.l(b,c) we show the vector field Vs k in the vicinity 
of the FP r>,i = (0.5098, -0.6084) and r Fj2 = (0.6216,0.6059) for the second iterate of the Ikeda 
map which represent a period two cycle of the map. In Fig. 1(b) 200 points in the neighbourhood 
of rp,i are plotted. The stabilization matrix Afc = 0.02 Ci was used. In Fig. 1(c) we use the same 
number of points as well as the same stabilization matrix Afe to illustrate the organization of the 
vector field around rp,2- Again the positions of the FP are indicated by crosses. From Figs.l(b,c) it is 
evident that the flow of the vector field Vs k can exhibit sharp, but still smooth, turns in the immediate 
neighbourhood of the FP. The size of these regions of nonlincarity in general decreases with increasing 
period of the cycle considered. 

The fact that the stabilized FP represent the centers of the flow of the vector field Vs k possesses 
a counterpart in the original chaotic system U which we shall discuss briefly in the following. If a 
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chaotic trajectory approaches a FP or more precisely enters the linear neighbourhood of the FP it 
can be shown that the maximum of the deflection of the trajectory in the linear regime occurs for the 
position closest to the FP. The deflection of the incident trajectory can be considered as a turning 
of the trajectory and the point of maximum deflection (or curvature for continuous systems) can 
be considered as a generalization of the notion of a turning point in one dimension (see refs. [^6[ ). 
Different chaotic trajectories approaching the FP from different directions in space are deflected to 
different directions of their outgoing manifold. The common feature of chaotic trajectories in the 
neighbourhood of a FP is the fact that the above-mentioned deflection or turning occurs with respect 
to any direction of space. If the chaotic system U is transformed via the corresponding stabilizing 
transformation Z^, the resulting dynamical system possesses centers of the corresponding flow of 
the transformed system at the positions of the FP. can therefore be considered as a transformation 
from spreading to focusing flow. 

The above-described method of detecting the unstable FP (periodic orbits) for a given chaotic 
dynamical system involves the parameter A which has to be chosen sufficiently small in order to 
transform the unstable FP of the original chaotic system to a stable one via the transformations L^. 
With increasing period of the FP to be located, the parameter A has to be chosen increasingly smaller 
in order to achieve stabilization. However we should not choose A too small since the convergence 
of the iterated transformed dynamical laws (see eq.(2)) to the FP will then become very slow and 
involves a waste of computer time. There is a simple way to avoid this tuning of the parameter A 
which makes our approach independent of the parameter A. Taking the limit A — > in eq.(2) we 
obtain 

(n+i - n) = ?=Ck (/ - (r) _ ^ (6) 

A — >o A 

This equation represents the continuous formulation of our transformed discrete dynamical systems in 
eq.(2). Its solutions possess the same nice properties as those of eq.(2) in the sense that any stabilized 
trajectory is running to a stable FP. The important difference however is that the solutions of eq.(6) do 
not depend on the parameter A and it can therefore easily be solved by using any suitable integration 
scheme like, for example, a Predictor-Corrector integrator. 

This means in particular that the basin of attraction for a certain FP is now independent of the 
choice of the parameter A and represents a characteristics of the originally unstable FP and the 
transformation Lk- In Figure 2 we show a typical ensemble of trajectories which are solutions to 
eq.(6) for a given set of starting points for the case of the stabilized period two orbit of the Ikeda map 
(see section 3). Part of the trajectories are running towards one point of the period two orbit and 
the other part towards the second point of the orbit. The boarder line of the two basins of attraction 
represents an unstable one dimensional manifold which is the only region of instability concerning the 
coordinate range of the Ikeda attractor. 

We emphasize that the above described method for the detection of the unstable periodic orbits can 
be applied to any chaotic dynamical system independent of its dimension and individual properties. 
The advantages of our approach are obvious: given any chaotic dynamical system with a skeleton 
of unstable periodic orbits we can transform these orbits into stable periodic orbits of a different 
transformed dissipative dynamical system. A large class of periodic orbits is thereby stabilized by 
the same linear transformation. The positions of the periodic orbits can then be obtained simply by 
iterating the transformed system or solving the differential equations of its continuous formulation. 
In this way the FP can be determined with in principle arbitrary accuracy. 
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Finally we mention that our method is by no means restricted to discrete dynamical systems. 
Periodic orbits in continuous dynamical systems can be detected by using the Poincare map which is 
again a discrete map representing the original continuous dynamical system in a chosen hyperspace. 
It is hereby not necessary for the Poincare map to be given in an analytical form. The transformed 
dynamical systems can be obtained from eqs.(2,6) for a numerically given Poincare map. All the 
quantities given in these equations are defined in the Poincare section which is a (N-l)-dimensional 
space and periodic orbits represent a finite number of points in this space i.e. are of integer period in 
the Poincare section. It is therefore not necessary to know the period of the orbits in the continuous 
system in order to detect them in the Poincare section. 

III. APPLICATION TO CHAOTIC DYNAMICAL SYSTEMS 

We applied our method to several 2D iterative maps like, for example, the Henon, 2D-logistic and 
in particular the more complicated Ikeda map. In the following we will discuss in some detail our 
procedure and the corresponding results for the Ikeda attractor. 

To demonstrate the reliability as well as efficiency of our method we use it to calculate the unstable 
periodic orbits of the Ikeda attractor. The underlying dynamical law is the Ikeda map |l5| used in 
nonlinear optics to describe the response of a 2-level homogeneous absorber in a ring cavity to a 
constant incident light wave. The map is given by: 

Ui : x n+ i = 1 + 0.9 (x n cos w n - y n sin w n ) 

y n+ i = 0.9 (x n sin w n + y n cos w n ) w n =0.4- — — ^— — - (7) 

The Ikeda attractor is fully chaotic and embedded in two dimensional coordinate space. Our con- 
siderations for the N — 2 case described in section 2 can therefore be applied directly. In detail we 
proceed as follows: Using the set of stabilization matrices Cj, i = 1, .., 5 given in section 2 we construct 
according to eqs.(2,6) the transformed systems which possess the corresponding stable FP. Starting 
from an initial point we next produce iteratively a chaotic trajectory on the attractor which provides 
a set of starting points for the trajectories of the transformed dynamical systems which are solutions 
to eqs.(2,6). If a certain trajectory converges to a point fp within a given accuracy then this point is 
recorded as a FP of the Ikeda map. We repeated this procedure for the higher iterates of the Ikeda 
map up to 13 — th iterate. The number of starting points on the attractor needed to obtain the 
periodic orbits of period 13 is roughly 5000. This is particularly impressive if we take into account 
that some of the period 13 orbits differ only from their fourth significant digit on and that the number 
of FP detected for period 13 is 2522. We are therefore in the position to resolve a large number of 
close lying periodic orbits. Again we emphasize that we have no rigorous proof for the completeness 
of the detected orbits but empirical evidence. 

The origin of the success of our stabilization transformation is its global character. Thus we need 
only a very coarse grained lattice of initial conditions to cover the attractor. For the discrete formu- 
lation of our transformed dynamical systems (see eq.(2)) the values of the parameter A necessary to 
achieve stabilization vary from 10 _1 for orbits of low period to a few times 10~ 5 for periodic orbits 
of period 13. Another important feature of our approach is the high accuracy which is controlled 
by the convergence of the corresponding trajectories. Our calculations yield a relative accuracy of 
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10 -13 . Having determined the coordinates of the cycles we use the Ikeda map in eq.(|7]) to derive the 
corresponding stability coefficients. 

The corresponding results are presented in table I. The first row shows the number of unstable 
cycles for a given period p and in the second row the total number of cycle points of order N is given. 
Starting with one period 2, two period 3 and three period 4 cycles their number increases strongly 
with increasing period and finally for period 12 and 13 we obtain 110 and 194 cycles, respectively. 
With the help of the obtained orbits and their stability properties we are in the position to determine 
characteristic quantities of the attractor related to its degree of chaoticity or dynamical and geometrical 
complexity. To get an impression of the quality of the covering (representation) of the Ikeda attractor 
through the unstable periodic orbits we illustrate in Fig. 3 (a) the set of all the FP of the Ikeda map and 
its higher iterates (up to the 13 — th iterate). For comparison we present in Fig. 3(b) a typical chaotic 
trajectory on the attractor which possesses the same total number of points (5627). Although the 
overall picture looks very similar a more careful comparison of Figs.3(a,b) reveals major differences in 
the local density of points. If we compare the distribution of the FP shown in Fig. 3(a) with the Ikeda 
attractor of Fig. 3(b) we realize that there exist regions within the support of the attractor which are 
not visited by periodic cycles up to period 13 at all. In addition there exist regions of high density 
of FP (see Fig. 3(a)) which do not correspond to regions of high density for the attractor itself (see 
Fig.3(b)). 

In order to perform a quantitative analysis of the Ikeda attractor in terms of the unstable cycles let 
us first calculate the measure of the exponential increase (with increasing period) of their number or 
in other words the so called topological entropy defined through |l7j : 

St = lim — lnn(p) (8) 

p — >oo p 

where nip) denotes the number of FP of Uf. The values of St are presented in table I (third row) 
and they can be used as a rough guide for the completeness of the number of periodic orbits found 
for a given period. An average measure of the strangeness of the attractor is given in terms of the 
average Lyapunov exponents defined through M : 

h e ,c = lim — —— -lniJ, e c (fjp) (9) 

where fj p denotes the j-th FP of the p-th iterate of the map U] and /i e (f\, p ) (/x c (f\, p )) is the expanding 
(contracting) eigenvalue of the stability matrix at this point. The two exponents correspond to an 
average expanding rate (h e ) and an average contracting rate (h c ), respectively. The results as a 
function of period p are presented in the fourth row of Table I. The strong fluctuations in successive 
terms of this expansion can be seen best in Fig. 4(a) where the two Lyapunov exponents are shown as 
a function of the period p. To characterize the strange attractor geometrically we follow [Q and try to 
cover the attractor with slabs of length 1 and width — tV^t . In this case the corresponding Hausdorff 
dimension D Q can be found by solving the equation: 

EMr^-^l (10) 

3 

The numerical results with increasing period p are given in the fifth and last row of table I. The 
corresponding graphical illustration is shown in Fig. 4(b). Again strong oscillations are present even 
for higher periods. 
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Summarizing our analysis of the Ikeda attractor we observe a satisfactory convergence of the topo- 
logical entropy up to period 13 but strongly oscillating values for the fractal dimension as well as 
average Lyapunov exponents. 

IV. CONCLUSIONS 

The main objective of the present paper is the development of a general method to detect the 
unstable periodic orbits in chaotic dynamical systems. The central idea is to convert the unstable 
fixed points (periodic orbits) to stable ones without changing their location in space. We have shown 
that this can be achieved by a set of linear transformations, i.e. transformation matrices, which allow 
the stabilization of any configuration of unstable hyperbolic fixed points for a given dynamical system. 
The numerical calculation of the stable fixed points in the transformed dissipative dynamical systems 
is done by either a simple iteration procedure or by solving the corresponding continuous version of 
the transformed dynamical laws. 

The above approach allows a straightforward application to any analytically or numerically given 
dynamical system. Apart from its general applicability the advantages of our method are the following. 
It is of universal character in the sense that no previous knowledges about the topological or dynamical 
behaviour, except the presence of chaos, are required. Our method is by no means restricted to lower 
dimensional (ID or 2D) systems but can in principle yield periodic orbits for dynamical systems 
of any dimension. It allows an efficient convergence to highly accurate values and requires at the 
same time only a small set of initial conditions to cover the strange attractor: a coarse grained 
covering is sufficient to detect the periodic orbits which are within our empirical procedure assumed 
to be complete. Nevertheless even close lying periodic orbits of higher periods can be resolved and, 
therefore, distinguished through the evolution of the stabilized system. The underlying stabilization 
transformation possesses an appealing geometrical interpretation. While the dynamics of the original 
chaotic system in the vinicity of the fixed points is characterized by a 'turning' of trajectories |l6f| 
in any direction of coordinate space the stabilized system belongs to a vector field which is centered 
around specified configurations of sinks/sources. The positions of those sinks/sources are identical 
with the positions of the fixed points. We remark that our method can also be applied to cases where 
the dynamical system is not fully chaotic but consists of a chaotic sea filled with islands of regularity. 

The above method is not restricted to discrete dynamical systems: periodic orbits in continuous 
dynamical systems can be detected by using the Poincare map which is again a discrete map repre- 
senting the original continuous dynamical system in a certain hyperspace. As a further perspective 
and application of the developed scheme we mention the possibility of solving nonlinear equations in 
general through the iteration of a suitably stabilized version of the original equations. Thereby we 
translate the problem of finding the roots of a set of nonlinear equations fix) = to that of finding 
the fixed points of the dynamical law F(x) — f(x) + x. Applying the above method to get a stabilized 
version of F we obtain then by simple iteration the fixed points of F which correspond to the roots 
off. 

Acknowledgments 

The European Community (F.K.D.) is gratefully acknowledged for financial support. 



10 



[1] D. Auerbach, P. Cvitanovic, J.-P. Eckmann, G. H. Gunaratne and I. Procaccia, Phys. Rev. Lett. 58, 2387 
(1987) 

[2] P. Cvitanovic, Phys. Rev. Lett. 61, 2729 (1988) 

[3] C. Grebogi, E. Ott, J. A. Yorke, Phys. Rev. A 37, 1711 (1988) 

[4] E. Ott, T. Sauer and J. A. Yorke, Phys. Rev. A 39, 4212 (1989) 

[5] D. P. Lathrop and A. J. Kostelich, Phys. Rev. A 40, 4028 (1989) 

[6] E. Ott, C. Grebogi and J. A. Yorke, Phys. Rev. Lett. 64, 1196 (1990). 

[7] M.C. Gutzwiller, " Chaos in Classical and Quantum Mechanics" , Springer- Verlag 1990 

[8] E. Ott, "Chaos in dynamical systems", Cambridge University Press (1992), p. 326 and references therein. 

[9] J.P. Eckmann and D. Ruelle, Rev.Mod.Phys.57, 617 (1985) 

[10] O. Biham and W. Wenzel, Phys. Rev. Lett. 63, 819 (1989); Phys. Rev. A 42, 4639 (1990) 

[11] P. Grassberger, H. Kantz and U. Moenig, J.Phys.A22, 5217 (1989) 

[12] K. Hansen, Phys. Rev. E 52, 2388 (1995). 

[13] P. Cvitanovic, G. H. Gunaratne and I. Procaccia, Phys. Rev. A 38, 1503 (1988) 

[14] P.Schmelcher and F.K.Diakonos, Phys.Rev.Lett.78, 4733 (1997) 

[15] K. Ikeda, Opt. Comm. 30, 257 (1979); S. M. Hammel, C. K. R. T. Jones and J. V. Moloney, J. Opt. Soc. 
Am. B 4, 552 (1985). 

[16] F.K. Diakonos and P. Schmelcher, CHAOS 7, 239 (1997); acc. f. publ. in Int. J. Bif. and Chaos. 

[17] A. Katok, Pub. Math. IHES 51, 377 (1980). 



Figure Captions 
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FIG. 1. (a) The vector field Vs k belonging to the Ikeda map on its attractor. (b) The corresponding vector 
field of the second iterate of the Ikeda map around the (stabilized) fixed point (0.5098,-0.6084) belonging to 
the period two cycle of the Ikeda map. (c) The vector field around the other fixed point (0.6216,0.6059) of 
the period two cycle. The positions of the fixed points are indicated by crosses. 200 points have been used for 
each subfigure. 

FIG. 2. Lines which represent an ensemble of trajectories of the transformed map in eq.(2) for the stabiliza- 
tion of the period 2 orbit of the Ikeda attractor. The circles indicate the positions of the points of the period 
2 orbit and the line connecting the crosses indicates the boarder line between the two basins of attraction. 

FIG. 3. (a) The set of all the fixed points of the Ikeda map and its higher iterates up to period 13. (b) A 
typical chaotic trajectory on the Ikeda attractor with the same number of points (5627) as in (a). 

FIG. 4. (a) The average Lyapunov exponents as a function of the period p. (b) The fractal dimension D 
of the Ikeda attractor as a function of the period p. See section 3. 



Table Captions 

Table I: The number of cycles with period N, the total number of cycle points of order N, the 
topological entropy, the Lyapunov exponents as well as the fractal dimension for N=l,...,13 (see 
section 3). 
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N 


1 


2 


3 


4 


5 


6 


7 


8 


9 


10 


11 


12 


13 


Number of cycles 
with period N 


1 


1 


2 


3 


4 


7 


10 


14 


26 


45 


76 


110 


194 


Number of cycle 
points order N 


1 


3 


7 


15 


21 


51 


71 


127 


241 


473 


837 


1383 


2523 


Topological 
Entropy 


0.0 


0.549 


0.648 


0.677 


0.608 


0.655 


0.609 


0.606 


0.609 


0.616 


0.612 


0.603 


0.603 


Lyapunov 
Exponents 


0.365 
-0.453 


0.513 
-0.648 


0.547 
-0.686 


0.735 
-0.942 


0.597 
-0.793 


0.884 
-1.144 


0.490 
-0.647 


1.071 
-1.493 


0.854 
-1.220 


0.957 
-1.313 


0.688 
-0.954 


1.046 
-1.388 


0.468 
-0.639 


Fractal 
Dimension 


1.00 


1.68 


1.76 


1.95 


1.83 


2.04 


1.78 


2.06 


1.95 


2.02 


1.81 


2.13 


1.76 



